
#### Replication script for paper: ####
#### How does Polycentric Engagement Relate to Countries' NDC Ambition and Mitigation Policy Effort? ####
#### by Paula Castro, Marlene Kammerer, Axel Michaelowa ####
#### published in Global Environmental Politics 24:3, 2024 ####

#### Questions and comments: paula.castro@zhaw.ch
#### Script run using R version 4.4.1


## Clean workspace, load required libraries and data ####
rm(list= ls())

# Check if the packages that we need are installed
want = c("dplyr", "ggplot2", "ggthemes", "ggpubr", "grid", "gridExtra", 
         "Hmisc", "interactions", "jtools", "MASS", "moments", "purrr", "readxl", 
         "stringr", "sjmisc", "sjPlot", "texreg", "tibble", "tidyr", "vtable")
have = want %in% rownames(installed.packages())

# Install the packages that we miss
if (any(!have)) {install.packages(want[!have])}

# Load the packages
lapply(want, library, character.only = TRUE)

# Remove the objects we created
rm(have, want)


## Load dataset
load(file = "final_dataset.RData")


#### 1. ADD SOME NECESSARY VARIABLES ####

# Logged variables
final_dataset$log_co2_cap_2014 <- log10(final_dataset$co2_cap_2014)
final_dataset$log_gdp_cap_2014 <- log10(final_dataset$gdp_cap_ppp_2014)

final_dataset$log_co2_cap_ndcyear <- log10(final_dataset$co2_cap_ndcyear)
final_dataset$log_gdp_cap_ndcyear <- log10(final_dataset$gdp_cap_ppp_ndcyear)

final_dataset$log_co2_cap_b4ndcyear <- log10(final_dataset$co2_cap_b4ndcyear)
final_dataset$log_gdp_cap_b4ndcyear <- log10(final_dataset$gdp_cap_ppp_b4ndcyear)

theme_set(theme_sjplot())


#### 2. MAIN RESULTS ####

##### Table 1: Variables and descriptive statistics ####

final_variables <- subset(final_dataset, select = c(cty_code, ambition_excl_uncond_cap_gdr, policy_density_2014, policy_density_strong_2014, instrument_density_strong_2014, 
                                                memberships_state, memberships_sub_nonstate, log_gdp_cap_b4ndcyear, fossil_rents_b4ndcyear, ndgain_b4ndcyear, br_dem_b4ndcyear, annex))

st(final_variables, labels = c("Ambition",
                               "Policy density",
                               "Strong policy density",
                               "Strong instrument density",
                               "State engagement",
                               "Nonstate engagement",
                               "GDP/cap (logged)",
                               "Fossil fuel rents",
                               "ND-Gain index",
                               "Democracy",
                               "Annex I membership"),
   out = "csv", file = "Table1_summary_stats.csv")


##### Table 2: Regression results for NDC ambition ####

ambition_postparis_1 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_2 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_3 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_4 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")

ambition_postparis_1lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_2lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_3lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_4lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(ambition_postparis_1, ambition_postparis_2, ambition_postparis_3, ambition_postparis_4),
        file = "Table2_Ambition_GLM_PostParis.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition" = 1:4))

htmlreg(list(ambition_postparis_1lm, ambition_postparis_2lm, ambition_postparis_3lm, ambition_postparis_4lm),
        file = "Table2_Ambition_LM_PostParis.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition" = 1:4))
# lm models: results are exactly the same; the output produces the adjusted R-squared.


##### Table 3: Regression results for mitigation policy effort ####

policydensity_postparis_1 <- glm.nb(policy_density_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_2 <- glm.nb(policy_density_2014 ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_3 <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_4 <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

policydensitystr_postparis_1 <- glm.nb(policy_density_strong_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_2 <- glm.nb(policy_density_strong_2014 ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_3 <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_4 <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

instrumentdensitystr_postparis_1 <- glm.nb(instrument_density_strong_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_2 <- glm.nb(instrument_density_strong_2014 ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_3 <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_4 <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(policydensity_postparis_1, policydensity_postparis_2, policydensity_postparis_3, policydensity_postparis_4, 
             policydensitystr_postparis_1, policydensitystr_postparis_2, policydensitystr_postparis_3, policydensitystr_postparis_4,
             instrumentdensitystr_postparis_1, instrumentdensitystr_postparis_2, instrumentdensitystr_postparis_3, instrumentdensitystr_postparis_4),
        file = "Table3_Policyeffort_GLM_PostParis.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Policy density" = 1:4, "Strong policy density" = 5:8, "Strong instrument density" = 9:12))


##### Figures 1-3: Interaction plots for mitigation policy effort ####

interact_plot(model = policydensity_postparis_2, pred = memberships_state, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure1a_interaction_policy_density_state_PostParis.png", width=12, height=9, units="cm")

interact_plot(model = policydensity_postparis_4, pred = memberships_sub_nonstate, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure1b_interaction_policy_density_nonstate_PostParis.png", width=12, height=9, units="cm")


interact_plot(model = policydensitystr_postparis_2, pred = memberships_state, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure2a_interaction_strong_policy_density_state_PostParis.png", width=12, height=9, units="cm")

interact_plot(model = policydensitystr_postparis_4, pred = memberships_sub_nonstate, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure2b_interaction_strong_policy_density_nonstate_PostParis.png", width=12, height=9, units="cm")


interact_plot(model = instrumentdensitystr_postparis_2, pred = memberships_state, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure3a_interaction_strong_instrument_density_state_PostParis.png", width=12, height=9, units="cm")

interact_plot(model = instrumentdensitystr_postparis_4, pred = memberships_sub_nonstate, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=11),
        axis.text = element_text(size=11),
        legend.title = element_text(size=11, face="bold"),
        legend.text = element_text(size=11),
        legend.position = "bottom")
ggsave(filename="Figure3b_interaction_strong_instrument_density_nonstate_PostParis.png", width=12, height=9, units="cm")


#### 3. APPENDIX TABLES AND GRAPHS ####

##### Appendix: 1. Descriptive statistics by group membership ####

final_variables <- subset(final_dataset, select = c(ambition_excl_uncond_cap_gdr, policy_density_2014, policy_density_strong_2014, instrument_density_strong_2014, 
                                                memberships_state, memberships_sub_nonstate, log_gdp_cap_b4ndcyear, fossil_rents_b4ndcyear, ndgain_b4ndcyear, br_dem_b4ndcyear, annex))

st(final_variables, labels = c("Ambition",
                               "Policy density",
                               "Strong policy density",
                               "Strong instrument density",
                               "State engagement",
                               "Nonstate engagement",
                               "GDP/cap (logged)",
                               "Fossil fuel rents",
                               "ND-Gain index",
                               "Democracy"), 
   group = 'annex', group.long = F, out = "csv", file = "Appendix1_summary_stats_bygroup.csv")


##### Appendix: 2. Correlation table ####

final_nona <- final_variables %>% drop_na() %>% dplyr::select(-annex)

correlation_matrix <- cor(final_nona)
write.csv(correlation_matrix, file = "Appendix2_correlation_matrix.csv")


##### Appendix: 3. Density plots of dependent variables by group membership ####

myPal <- get_colors("CUD Bright", 2) # tells me the colors that were used in the interaction plots above

# Ambition
ggplot(final_dataset, aes(x=ambition_excl_uncond_cap_gdr, fill=annex, color=annex)) +
  geom_density(alpha=.3) + scale_fill_manual(values = c("#49b7fc", "#ff7b00")) + scale_color_manual(values = c("#49b7fc", "#ff7b00")) +
  labs(fill = "UNFCCC group membership") +
  labs(color = "UNFCCC group membership") +
  xlab("NDC ambition") +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom",
        line = element_line(linewidth = .7))

ggsave(filename = "Appendix3_Density_plot_NDC_ambition.pdf", width=12, height=9, units="cm")


# Policy density
ggplot(final_dataset, aes(x=policy_density_2014, fill=annex, color=annex)) +
  geom_density(alpha=.3) + scale_fill_manual(values = c("#49b7fc", "#ff7b00")) + scale_color_manual(values = c("#49b7fc", "#ff7b00")) +
  labs(fill = "UNFCCC group membership") +
  labs(color = "UNFCCC group membership") +
  xlab("Policy density") +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom",
        line = element_line(linewidth = .7))

ggsave(filename = "Appendix3_Density_plot_policy_density.pdf", width=12, height=9, units="cm")


# Strong policy density
ggplot(final_dataset, aes(x=policy_density_strong_2014, fill=annex, color=annex)) +
  geom_density(alpha=.3) + scale_fill_manual(values = c("#49b7fc", "#ff7b00")) + scale_color_manual(values = c("#49b7fc", "#ff7b00")) +
  labs(fill = "UNFCCC group membership") +
  labs(color = "UNFCCC group membership") +
  xlab("Strong policy density") +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom",
        line = element_line(linewidth = .7))

ggsave(filename = "Appendix3_Density_plot_strong_policy_density.pdf", width=12, height=9, units="cm")


# Strong instrument density
ggplot(final_dataset, aes(x=instrument_density_strong_2014, fill=annex, color=annex)) +
  geom_density(alpha=.3) + scale_fill_manual(values = c("#49b7fc", "#ff7b00")) + scale_color_manual(values = c("#49b7fc", "#ff7b00")) +
  labs(fill = "UNFCCC group membership") +
  labs(color = "UNFCCC group membership") +
  xlab("Strong instrument density") +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom",
        line = element_line(linewidth = .7))

ggsave(filename = "Appendix3_Density_plot_strong_instrument_density.pdf", width=12, height=9, units="cm")


##### Appendix: 4. Policy types in the CCLW dataset ####

# The data on policies and instruments used for this publication was obtained from the ‘Climate Change Laws of the World’ 
# interface for the Climate Policy Radar Database, https://climate-laws.org and made available under the 
# Creative Commons CC-BY licence.”

# Citation:
# Grantham Research Institute at the London School of Economics and Climate Policy Radar (2023): 
# Climate Change Laws of the World. https://climate-laws.org and https://app.climatepolicyradar.org/search
# Data downloaded on 26.08.2023. 

# The data was pre-processed to:
# Gather information about the latest change in a policy (this can be, for example, entry into force, an amendment, or issuance of an additional policy regulating aspects of a law),
# so that only those policies adopted or modified from 2014 onward are kept
# Keep only the main observation for each CPR policy family; 
# e.g., if there is a law and an ordinance regulating that law, both are separate observations that belong to the same policy family as per the CCLW dataset; 
# in this case, only the observation relating to the law is kept
# Keep only those policies that address mitigation (even though they may also address adaptation or other aspects of climate policy)

load(file = "CCLW_mitigation_policies.RData")

# Classify policies into strong (policies that actually implement soemthing) and weak ones (policies aimed at planning future actions, but that do not yet implement them)
CCLW_mitigation_policies$type_broad <- ifelse(str_detect(CCLW_mitigation_policies$Document_Type, "Act|Decision|Decree|Directive|EU Decision|EU Directive|EU Regulation|Executive Order|Law|Order|Ordinance|Regulation|Resolution|Royal Decree|Rules"), "Strong", "Weak")
CCLW_mitigation_policies$type_broad[CCLW_mitigation_policies$Document_Type == "Action Plan"] <- "Weak"

# Table of strong policies
table(CCLW_mitigation_policies$Document_Type[CCLW_mitigation_policies$type_broad == "Strong"]) %>% 
  as.data.frame() %>% 
  arrange(desc(Freq))

# Table of weak policies
table(CCLW_mitigation_policies$Document_Type[CCLW_mitigation_policies$type_broad == "Weak"]) %>% 
  as.data.frame() %>% 
  arrange(desc(Freq))


##### Appendix: 5. Instrument types in the CCLW dataset ####

# Split policies dataset into individual policy instruments 
mitigation_policies_splitted <- CCLW_mitigation_policies %>%
  transform(Instruments = strsplit(Instruments, ";", fixed = TRUE)) %>% # splits the column "Instruments", creating a list of all instruments within the column
  unnest(Instruments) # unnests the list into different rows
# from 2539 to 4885 observations

# Classify policy instruments into strong and weak ones
mitigation_policies_splitted$instrument_type <- ifelse(str_detect(mitigation_policies_splitted$Instruments, "Economic|Regulation|Direct Investment"), "Strong", "Weak")

# Get rid of early warning systems, as not mitigation-related
mitigation_policies_splitted <- subset(mitigation_policies_splitted, Instruments != "Early warning systems|Direct Investment")
# now 4858 policy instruments

# Table of strong policy instruments
table(mitigation_policies_splitted$Instruments[mitigation_policies_splitted$instrument_type == "Strong"]) %>% 
  as.data.frame() %>% 
  arrange(desc(Freq))

# Table of weak policy instruments
table(mitigation_policies_splitted$Instruments[mitigation_policies_splitted$instrument_type == "Weak"]) %>% 
  as.data.frame() %>% 
  arrange(desc(Freq))


##### Appendix: 6. Additional regression results with alternative measures of NDC ambition ####

## A. Per capita approach
ambition_postparis_1_cap <- glm(ambition_excl_uncond_cap_cap ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_2_cap <- glm(ambition_excl_uncond_cap_cap ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_3_cap <- glm(ambition_excl_uncond_cap_cap ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_4_cap <- glm(ambition_excl_uncond_cap_cap ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")

ambition_postparis_1lm_cap <- lm(ambition_excl_uncond_cap_cap ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_2lm_cap <- lm(ambition_excl_uncond_cap_cap ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_3lm_cap <- lm(ambition_excl_uncond_cap_cap ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_4lm_cap <- lm(ambition_excl_uncond_cap_cap ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(ambition_postparis_1_cap, ambition_postparis_2_cap, ambition_postparis_3_cap, ambition_postparis_4_cap),
        file = "Appendix6a_Ambition_GLM_PostParis_capita.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition, per capita approach" = 1:4))

htmlreg(list(ambition_postparis_1lm_cap, ambition_postparis_2lm_cap, ambition_postparis_3lm_cap, ambition_postparis_4lm_cap),
        file = "Appendix6a_Ambition_LM_PostParis_capita.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE, custom.header = list("Ambition, per capita approach" = 1:4))


## B. Cumulative per capita approach
ambition_postparis_1_cumcap <- glm(ambition_excl_uncond_cap_cumcap ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_2_cumcap <- glm(ambition_excl_uncond_cap_cumcap ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_3_cumcap <- glm(ambition_excl_uncond_cap_cumcap ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_4_cumcap <- glm(ambition_excl_uncond_cap_cumcap ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")

ambition_postparis_1lm_cumcap <- lm(ambition_excl_uncond_cap_cumcap ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_2lm_cumcap <- lm(ambition_excl_uncond_cap_cumcap ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_3lm_cumcap <- lm(ambition_excl_uncond_cap_cumcap ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_4lm_cumcap <- lm(ambition_excl_uncond_cap_cumcap ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(ambition_postparis_1_cumcap, ambition_postparis_2_cumcap, ambition_postparis_3_cumcap, ambition_postparis_4_cumcap),
        file = "Appendix6b_Ambition_GLM_PostParis_cumulative_capita.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition, cumulative per capita" = 1:4))

htmlreg(list(ambition_postparis_1lm_cumcap, ambition_postparis_2lm_cumcap, ambition_postparis_3lm_cumcap, ambition_postparis_4lm_cumcap),
        file = "Appendix6b_Ambition_LM_PostParis_cumulative_capita.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE, custom.header = list("Ambition, cumulative per capita" = 1:4))


## C. Capability approach
ambition_postparis_1_capab <- glm(ambition_excl_uncond_cap_capab ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_2_capab <- glm(ambition_excl_uncond_cap_capab ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_3_capab <- glm(ambition_excl_uncond_cap_capab ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")
ambition_postparis_4_capab <- glm(ambition_excl_uncond_cap_capab ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset, family = "gaussian")

ambition_postparis_1lm_capab <- lm(ambition_excl_uncond_cap_capab ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_2lm_capab <- lm(ambition_excl_uncond_cap_capab ~ memberships_state * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_3lm_capab <- lm(ambition_excl_uncond_cap_capab ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
ambition_postparis_4lm_capab <- lm(ambition_excl_uncond_cap_capab ~ memberships_sub_nonstate * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(ambition_postparis_1_capab, ambition_postparis_2_capab, ambition_postparis_3_capab, ambition_postparis_4_capab),
        file = "Appendix6c_Ambition_GLM_PostParis_capabilities.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition, capabilities approach" = 1:4))

htmlreg(list(ambition_postparis_1lm_capab, ambition_postparis_2lm_capab, ambition_postparis_3lm_capab, ambition_postparis_4lm_capab),
        file = "Appendix6c_Ambition_LM_PostParis_capabilities.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE, custom.header = list("Ambition, capabilities approach" = 1:4))


##### Appendix: 7. Additional regression results with 3 country categories, including LDCs ####

policydensity_postparis_1 <- glm.nb(policy_density_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_2_ldcs <- glm.nb(policy_density_2014 ~ memberships_state * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_3 <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_postparis_4_ldcs <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

policydensitystr_postparis_1 <- glm.nb(policy_density_strong_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_2_ldcs <- glm.nb(policy_density_strong_2014 ~ memberships_state * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_3 <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_postparis_4_ldcs <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

instrumentdensitystr_postparis_1 <- glm.nb(instrument_density_strong_2014 ~ memberships_state + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_2_ldcs <- glm.nb(instrument_density_strong_2014 ~ memberships_state * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_3 <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_postparis_4_ldcs <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate * ldc_cats + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

screenreg(list(policydensity_postparis_1, policydensity_postparis_2_ldcs, policydensity_postparis_3, policydensity_postparis_4_ldcs,
               policydensitystr_postparis_1, policydensitystr_postparis_2_ldcs, policydensitystr_postparis_3, policydensitystr_postparis_4_ldcs,
               instrumentdensitystr_postparis_1, instrumentdensitystr_postparis_2_ldcs, instrumentdensitystr_postparis_3, instrumentdensitystr_postparis_4_ldcs))

interact_plot(model = policydensity_postparis_2_ldcs, pred = memberships_state, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Policy density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7a_Interaction_policy_density_state_PostParis_ldcs.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensity_postparis_4_ldcs, pred = memberships_sub_nonstate, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Policy density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7a_Interaction_policy_density_nonstate_PostParis_ldcs.pdf", width=12, height=9, units="cm")


interact_plot(model = policydensitystr_postparis_2_ldcs, pred = memberships_state, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong policy density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7b_Interaction_strong_policy_density_state_PostParis_ldcs.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensitystr_postparis_4_ldcs, pred = memberships_sub_nonstate, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong policy density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7b_Interaction_strong_policy_density_nonstate_PostParis_ldcs.pdf", width=12, height=9, units="cm")


interact_plot(model = instrumentdensitystr_postparis_2_ldcs, pred = memberships_state, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong instrument density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7c_Interaction_strong_instrument_density_state_PostParis_ldcs.pdf", width=12, height=9, units="cm")

interact_plot(model = instrumentdensitystr_postparis_4_ldcs, pred = memberships_sub_nonstate, modx = ldc_cats, interval = TRUE,
              colors = c("deepskyblue", "grey50", "darkorange1"),
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong instrument density", legend.main = "Group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix7c_Interaction_strong_instrument_density_nonstate_PostParis_ldcs.pdf", width=12, height=9, units="cm")


##### Appendix: 8. Additional regression results with explanatory variables (polycentric engagement) measured before December 2020 (before latest NDC submission) ####

data_latestNDCs <- subset(final_dataset, latest_ndc_year > 2019)
data_latestNDCs <- subset(data_latestNDCs, cty_code %nin% c("MNG", "CUB", "JAM", "RWA", "CHL", "MDA")) # these countries submitted their latest NDC in 2020, but before December
# 160 observations

ambition_pre2020_1 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_state_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs, family = "gaussian")
ambition_pre2020_2 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_state_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs, family = "gaussian")
ambition_pre2020_3 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs, family = "gaussian")
ambition_pre2020_4 <- glm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs, family = "gaussian")

ambition_pre2020_1lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_state_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs)
ambition_pre2020_2lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_state_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs)
ambition_pre2020_3lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs)
ambition_pre2020_4lm <- lm(ambition_excl_uncond_cap_gdr ~ memberships_sub_nonstate_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = data_latestNDCs)

policydensity_pre2020_1 <- glm.nb(policy_density_2014 ~ memberships_state_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_pre2020_2 <- glm.nb(policy_density_2014 ~ memberships_state_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_pre2020_3 <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensity_pre2020_4 <- glm.nb(policy_density_2014 ~ memberships_sub_nonstate_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

policydensitystr_pre2020_1 <- glm.nb(policy_density_strong_2014 ~ memberships_state_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_pre2020_2 <- glm.nb(policy_density_strong_2014 ~ memberships_state_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_pre2020_3 <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
policydensitystr_pre2020_4 <- glm.nb(policy_density_strong_2014 ~ memberships_sub_nonstate_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

instrumentdensitystr_pre2020_1 <- glm.nb(instrument_density_strong_2014 ~ memberships_state_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_pre2020_2 <- glm.nb(instrument_density_strong_2014 ~ memberships_state_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_pre2020_3 <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate_2020 + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)
instrumentdensitystr_pre2020_4 <- glm.nb(instrument_density_strong_2014 ~ memberships_sub_nonstate_2020 * annex + log_gdp_cap_b4ndcyear + fossil_rents_b4ndcyear + ndgain_b4ndcyear + br_dem_b4ndcyear, data = final_dataset)

htmlreg(list(ambition_pre2020_1lm, ambition_pre2020_2lm, ambition_pre2020_3lm, ambition_pre2020_4lm),
        file = "Appendix8_Ambition_LM_Pre2020.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition" = 1:4))

htmlreg(list(ambition_pre2020_1, ambition_pre2020_2, ambition_pre2020_3, ambition_pre2020_4,
             policydensity_pre2020_1, policydensity_pre2020_2, policydensity_pre2020_3, policydensity_pre2020_4,
             policydensitystr_pre2020_1, policydensitystr_pre2020_2, policydensitystr_pre2020_3, policydensitystr_pre2020_4, 
             instrumentdensitystr_pre2020_1, instrumentdensitystr_pre2020_2, instrumentdensitystr_pre2020_3, instrumentdensitystr_pre2020_4),
        file = "Appendix8_All_Pre2020.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition" = 1:4, "Policy density" = 5:8, "Strong policy density" = 9:12, "Strong instrument density" = 13:16))


interact_plot(model = policydensity_pre2020_2, pred = memberships_state_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8a_Interaction_policy_density_state_Pre2020.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensity_pre2020_4, pred = memberships_sub_nonstate_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8a_Interaction_policy_density_nonstate_Pre2020.pdf", width=12, height=9, units="cm")


interact_plot(model = policydensitystr_pre2020_2, pred = memberships_state_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8b_Interaction_strong_policy_density_state_Pre2020.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensitystr_pre2020_4, pred = memberships_sub_nonstate_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8b_Interaction_strong_policy_density_nonstate_Pre2020.pdf", width=12, height=9, units="cm")


interact_plot(model = instrumentdensitystr_pre2020_2, pred = memberships_state_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8c_Interaction_strong_instrument_density_state_Pre2020.pdf", width=12, height=9, units="cm")

interact_plot(model = instrumentdensitystr_pre2020_4, pred = memberships_sub_nonstate_2020, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix8c_Interaction_strong_instrument_density_nonstate_Pre2020.pdf", width=12, height=9, units="cm")


##### Appendix: 9. Additional regression results with explanatory variables (polycentric engagement) in the Pre-Paris period ####

ambition_preparis_1 <- glm(ambition_excl_uncond_cap_gdr ~ climate_memberships_rowan + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset, family = "gaussian")
ambition_preparis_2 <- glm(ambition_excl_uncond_cap_gdr ~ climate_memberships_rowan * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset, family = "gaussian")
ambition_preparis_3 <- glm(ambition_excl_uncond_cap_gdr ~ tcgs + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset, family = "gaussian")
ambition_preparis_4 <- glm(ambition_excl_uncond_cap_gdr ~ tcgs * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset, family = "gaussian")

ambition_preparis_1lm <- lm(ambition_excl_uncond_cap_gdr ~ climate_memberships_rowan + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
ambition_preparis_2lm <- lm(ambition_excl_uncond_cap_gdr ~ climate_memberships_rowan * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
ambition_preparis_3lm <- lm(ambition_excl_uncond_cap_gdr ~ tcgs + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
ambition_preparis_4lm <- lm(ambition_excl_uncond_cap_gdr ~ tcgs * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)

policydensity_preparis_1 <- glm.nb(policy_density_2014 ~ climate_memberships_rowan + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensity_preparis_2 <- glm.nb(policy_density_2014 ~ climate_memberships_rowan * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensity_preparis_3 <- glm.nb(policy_density_2014 ~ tcgs + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensity_preparis_4 <- glm.nb(policy_density_2014 ~ tcgs * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)

policydensitystr_preparis_1 <- glm.nb(policy_density_strong_2014 ~ climate_memberships_rowan + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensitystr_preparis_2 <- glm.nb(policy_density_strong_2014 ~ climate_memberships_rowan * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensitystr_preparis_3 <- glm.nb(policy_density_strong_2014 ~ tcgs + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
policydensitystr_preparis_4 <- glm.nb(policy_density_strong_2014 ~ tcgs * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)

instrumentdensitystr_preparis_1 <- glm.nb(instrument_density_strong_2014 ~ climate_memberships_rowan + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
instrumentdensitystr_preparis_2 <- glm.nb(instrument_density_strong_2014 ~ climate_memberships_rowan * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
instrumentdensitystr_preparis_3 <- glm.nb(instrument_density_strong_2014 ~ tcgs + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)
instrumentdensitystr_preparis_4 <- glm.nb(instrument_density_strong_2014 ~ tcgs * annex + log_gdp_cap_2014 + fossil_rents_2014 + ndgain_2014 + br_dem_2014, data = final_dataset)

htmlreg(list(ambition_preparis_1lm, ambition_preparis_2lm, ambition_preparis_3lm, ambition_preparis_4lm),
        file = "Appendix9_Ambition_LM_PreParis.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE)

htmlreg(list(ambition_preparis_1, ambition_preparis_2, ambition_preparis_3, ambition_preparis_4,
             policydensity_preparis_1, policydensity_preparis_2, policydensity_preparis_3, policydensity_preparis_4,
             policydensitystr_preparis_1, policydensitystr_preparis_2, policydensitystr_preparis_3, policydensitystr_preparis_4,
             instrumentdensitystr_preparis_1, instrumentdensitystr_preparis_2, instrumentdensitystr_preparis_3, instrumentdensitystr_preparis_4),
        file = "Appendix9_All_GLM_PreParis.doc", inline.css = FALSE, digits = 3, stars = c(0.01, 0.05, 0.1),
        doctype = TRUE, html.tag = TRUE, head.tag = TRUE, body.tag = TRUE, single.row = TRUE,
        custom.header = list("Ambition" = 1:4, "Policy density" = 5:8, "Strong policy density" = 9:12, "Strong instrument density" = 13:16))


interact_plot(model = policydensity_preparis_2, pred = climate_memberships_rowan, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9a_interaction_policy_density_state_PreParis.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensity_preparis_4, pred = tcgs, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9a_interaction_policy_density_nonstate_PreParis.pdf", width=12, height=9, units="cm")


interact_plot(model = policydensitystr_preparis_2, pred = climate_memberships_rowan, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9b_interaction_strong_policy_density_state_PreParis.pdf", width=12, height=9, units="cm")

interact_plot(model = policydensitystr_preparis_4, pred = tcgs, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong policy density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9b_interaction_strong_policy_density_nonstate_PreParis.pdf", width=12, height=9, units="cm")


interact_plot(model = instrumentdensitystr_preparis_2, pred = climate_memberships_rowan, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "State engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9c_interaction_strong_instrument_density_state_PreParis.pdf", width=12, height=9, units="cm")

interact_plot(model = instrumentdensitystr_preparis_4, pred = tcgs, modx = annex, interval = TRUE,
              int.type = "confidence", int.width = 0.95, outcome.scale = "link", x.label = "Nonstate engagement", y.label = "Strong instrument density", legend.main = "UNFCCC group membership", line.thickness = .7) +
  theme_bw() + 
  theme(axis.title = element_text(size=9),
        legend.title = element_text(size=9, face="bold"),
        legend.position = "bottom")
ggsave(filename="Appendix9c_interaction_strong_instrument_density_nonstate_PreParis.pdf", width=12, height=9, units="cm")
